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Abstract 

Fractional moments have been investigated by many authors to represent the 
density of univariate and bivariate random variables in different contexts. Frac- 
tional moments are indeed important when the density of the random variable has 
inverse power-law tails and, consequently, it lacks integer order moments. In this 
paper, starting from the Mellin transform of the characteristic function and by frac- 
tional calculus method we present a new perspective on the statistics of random 
variables. Introducing the class of complex moments, that include both integer and 
fractional moments, we show that every random variable can be represented within 
this approach, even if its integer moments diverge. Applications to the statistical 
characterization of raw data and in the representation of both random variables 
and vectors are provided, showing that the good numerical convergence makes the 
proposed approach a good and reliable tool also for practical data analysis. 

1 Introduction 

Fractional moments are very useful in dealing with random variables with power-law 
distributions, F(x) ~ \x\~ 11 , fi > 0, where F(x) is the distribution function. Indeed, in 
such cases, moments exist only if q < fi and integer order moments greater than 
/j diverge. Distributions of this type are encountered in a wide variety of contexts, see 
the extensive literature in JT) and ||2j where power-law statistics appear in the frame- 
work of anomalous diffusion in many fields of applied science. To name but a few 
examples, we mention the travel length distribution in human motion patterns and ani- 
mal search processes J3]|4j[5][6), fluctuations in plasma devices [7, 8|, scaling laws in 
polymer physics and its applications to gene regulation models |9 10 1, as well as the 
distribution of time scales in processes such as the motion of charge carriers in amor- 
phous semiconductors IfTTl . tracer dispersion in groundwater Ifl2l . or sticking times in 
turbulent flows [131. 
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The lack of moments is of course a great limitation for the characterization of such 
distributions, because many methods of statistical analysis fail, like for example, the 
representation of the characteristic function or the log-characteristic function by mo- 
ments or by cumulants, respectively. Then, it is necessary to resort to other distribution 
properties and a good option is given by so-called fractional moments which have been 
extensively studied in the seventies and recently attracted new interest both in theo- 
retical irT4l [T5l [TBI and in experimental settings [17|. In particular, in lfT31 and lfl6l 
the moment problem originally stated in terms of integer moments, has been extended 
to fractional moments where it has been shown that the knowledge of some fractional 
moments {X q ) improves significantly the convergence speed of the maximum entropy 
method for non-negative densities i.e. defined in [0, oo[, if a proper optimization pro- 
cedure is used in selecting the order q. This can be explained by the non-local nature 
of fractional moments that will be better highlighted in the following of the paper. 

Fractional moments of a non-negative random variable are expressible by Mellin 
transform of density and this fact has been widely used in literature principally in the 
field of the algebra of random variables. That is, Mellin transform is the principal 
mathematical tool to handle problems involving products and quotients of independent 
random variables. As the approach presented in this paper is based on Mellin transform, 
a brief outline of previous results is given in the following. 

The first author who developed a systematic method to express the density of the 
product of independent random variables by Mellin transform was Epstein. As con- 
sequence, in his pioneering work Ifl8l he recognizes that any density function can be 
obtained as Mellin convolution of two independent densities. It must be remarked that 
such a result was already obtained by lTT9l outside the framework of integral trans- 
forms. Many other studies on the use of Mellin transform extending Epstein's work 
have been reported providing explicit analytic forms for products of independent ran- 
dom variables with assigned densities. In 11201 examples on product and quotient of 
independent random variables with Rayleigh distribution and moments of Rice distri- 
bution are obtained by Mellin transform tables. An extensive use of such concepts and 
the state-of-art on the algebra of random variable may be found in the the book of CD . 

Another research direction on Mellin applications in probability is represented by 
the use of special functions like the H-Fox and the Maijer's G-functions, due principally 
to Mathai and co-workers. Such functions are indeed representable as Mellin-Barnes 
integrals of the product of gamma functions and are therefore suited to represent statis- 
tics of products and quotients of independent random variables whose fractional mo- 
ments are expressible as gamma or gamma related functions. Applications of special 
functions to statistics and probability theory can be found in 1221 23] [24] [25] [26] [27) 
and references therein. 

Whether fractional moments exists and how they are correlated to local properties 
of the characteristic function has been investigated in 11281 [29l l30l [3T1 . In particular, 
in [28 1 a relation is derived between the absolute moment of real order and integrals 
involving the Marchaud fractional derivative of the characteristic function. Other rel- 
evant studies on the existence of fractional (absolute) moments can also be found in 
G2] [33] [34). 

All the cited work share the common issue of using the Mellin transform because 
it naturally coincides with moments of the type (X 7 ~ l ) if the random variable X has a 
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density defined only in the positive domain. In this paper, following the previous results 
in [35 1, it is shown that applying the Mellin transform to the characteristic function, 
and not to the density, a sound representation of the statistics of general multivariate 
random variables is possible. In this way, exploiting the Hermitian nature of the char- 
acteristic function, a unique representation valid for every random variable and vectors, 
defined in the whole real range is presented. Further, by means of simple properties 
of fractional operators, a generalization of the well-known link between moments and 
the derivative of the characteristic function evaluated in zero is provided. Such a link is 
useful because it allows one to find existence criteria for fractional moments that can be 
derived in a simpler way with respect to [28 1 and clarifies the non-local nature of frac- 
tional moments. Moreover, this approach that starts from the characteristic function 
plainly applies to random variables defined in the whole domain without the need of 
partitioning. Such theoretical aspects go along with a better numerical treatment of the 
inverse Mellin transform involved in this approach as it will stressed in the following. 

In the application section, it will be indeed shown that the proposed approach is 
extremely useful in the statistical characterization of raw data and in the representation 
of both random variables and vectors with few computational effort. Thus, the good 
numerical convergence makes the approach proposed a good and reliable tool also for 
data analysis. We want to stress that the representation proposed in this paper can 
be applied to any Fourier pair functions used in statistical analysis and, very recently, 
applications to the representation of stochastic processes have been proposed in ||36l 

E2. 



2 Probabilistic Characterization of Random Variables 

A brief outline of the basic tools of probability will serve to frame the problem in hand 
and to introduce symbols. Let (S, M, f) be a probability space, where 3i is a collection 
of events defining a cr-algebra of subsets on the sample space S, and f a probability 
measure that assigns a number to each event on J[ and X : (S, — > (R, T>) a random 
variable, where S is a Borel set. F (x) = Pr (X < x) is the cumulative distribution and 
p (x) = dF/dx is the probability density function (PDF). The Fourier transform of p{x), 
denoted as Tp = <p(u), is the characteristic function (CF), that is 



f 

tp{u) — (exp(iuX)} — I exp (iux) p (x) dx (1) 

yj —OO 

with u e R, i = V-T and (■) means expectation, i.e. (g(X)} = g (x) p (x) dx. 
Expectations of the functions g (X) = X' with j - 1,2, provided they exist, give the 
integer moments of X, indicated by (X-i). These integer moments are related to the CF 
by the Taylor expansion 

OO j 

¥>(«) = £<(«)'>- (2) 

due to the property 

<(ixy>= dj(p(u)/du j \, (3) 
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Whilst eq.© might be useful to evaluate the CF from moments, its Fourier trans- 
form that expresses the relation between density and moments is not useful from the 
computational point of view, as it reads 



7=0 



where is the Dirac's delta. Then, it is not possible to use integer moments to 
represent the PDF by eq.©. It must be stressed that eq.© does not hold for every 
distribution. This is the case for distributions with heavy tails, and consequently with 
divergent integer moments. Further, the use of eq.(f2| is not always feasible, although 
appealing, because truncation of the Taylor expansion involved produces a divergent 
behavior of the CF at u — > oo. 

In a recent paper [35), it has been proposed to use a new class of moments, namely 
complex moments, to represent both the CF and the PDF of a random variable. The 
fundamentals of the generalization of eq.(f2]), eqs.© and (01 are briefly summarized 
hereinafter and it will be shown that such representation does not entail the drawback 
of eqs.(f2]l and eq.©. 

To this aim, let y — p + irj, with p > and rj e R, and let (x) and (D y ± f) 00 

denote the Riemann-Liouville (RL) fractional integral and derivative, while (p r ± f\ (x) 
is the Marchaud fractional derivative, namely 

(llf)(xA f ^ ^ e-\f(x + 0d{ (5a) 

( D ^)« d =nT^)sr rT/(Jc ^ )d ^ (5b) 

(Dlf) Qc/^—L-J Z-y- l (f(x + 0-f(x))df (5c) 

being F(y) the gamma function, and n the integer part of p+l. 

Fourier transform of fractional derivatives, in case < p < 1, (OH, p. 137) gives 

r(iif) = (+iurrf (6a) 
r (Dlf) = r (Dlf) = (Tiuy rf (6b) 

while inverse Fourier transform gives 

T- x (l y ± f) = (±xx)-VT- l f (7a) 

T- 1 (Dlf) = T- 1 (Dlf) = (±ixy r~ l f (7b) 
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For functions that are Fourier pairs, the previous equations are very useful for cal- 
culating the fractional operators in an easier way with respect to definitions (0. Indeed, 
by Fourier transform of eq.(0 



Xoo 
(+ixy y exp(iux)f(x)dx 
CO 

(D y ± Tf) (u) = (D y ± Tf) (u) = I (±ixY exp (iux) f (x) dx 



(8a) 



(8b) 



one obtains integrals that are not convolutions. Choosing as Fourier pair p(x) and (p{u), 
defined in eq.© and letting u = 0, it is easy to obtain the generalized form of eq.© 
from the latter, that can be written as 



(±i*rM*)d* = <(±uo- r > 

oo 

/CO 
(±bcy P (x)dx=((±ix)?) 
oo 



(9a) 



(9b) 



The set of complex moments is a natural generalization of integer moments as like 
as RL fractional differential operators generalize the classical differential calculus and, 
coherently eq.d9bl) coincides with eq.([3j, when y assumes integer values. The latter 
equations relate the behavior of the characteristic function and the existence of mo- 
ments in a more direct way than those proposed up to now in literature 1 28 1 . Moreover, 
in contrast to the local nature of the derivative and consequently of the integer moments, 
the complex moments are non-local as like the fractional derivatives. 

Whilst eqs.((9a| and d9bb descend from a proper application of the Fourier proper- 
ties of the RL operators, the generalization of eq.(f2| is not trivial; indeed, it relies on a 
generalized form of Taylor series. In ll35ll various generalizations of the Taylor expan- 
sion presented in literature have been analyzed, showing that the integral Taylor form 
based on inverse Mellin transform, sketched in the book of Samko et al. ( Il38ll . pp. 144- 
5), is the most useful to our scope. The starting point to find such a generalized Taylor 
form is to interpret the fractional derivative and integral defined in eq.© as a Mellin 
transform. In what follows we develop our proof just for a characteristic function <p(s) 
for the sake of simplicity, but the result obtained here is more general. Moreover, from 
this point on, we will use only the definition of the Marchaud fractional derivative be- 
cause same concepts follow plainly for the RL definition. Then, letting x = 0, eqs.® 
can be written as 

r(y)(f y ± <p)(0) = M{<p(+Z),y} (10a) 

T (-y) (Dl<p) (0) - M{<p (T0 - <p (0) , -y) (10b) 

where vVf{-,y} is the Mellin transform. Applying the inverse Mellin transform and 
recalling the condition i/?(0) = 1, two representations of the characteristic function 

^ /-»p+ioo 

<p(+u) = — T(y)(ll^)(Q) U -My (11a) 

£W J p-ioo 
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i r>p+ioo 

<^(+ M ) = l + — r(-y)(Di<p)(Q)urdy (lib) 

jp-ioo 

hold true, with u > and where integrals are performed along the imaginary axis with 
fixed real p that belongs to the so called fundamental strip of the Mellin transform of 
the function <p (u). 

The latter equations represent integral forms of the Taylor expansion for the CF, 
in the sense that from knowledge of the function (l r ± <p) (0) or (D y ± <p\ (0) it is possible 
to evaluate the function tp (u) by integration along the imaginary line with real part p 
(08). pp. 144-5). 

Such a representation gains an appealing touch once eqs.© are considered, leading 

to 

^ r*p+ioo 

^(+m)=^-J T(y){(±iX)~ y )u y dy (12a) 

i pp+ioo 

V (+ U ) = l + — r(- 7 )<(±iXF> M ?dy (12b) 

27T1 J p - ioo 

with u > 0. The latter equations are the integral extensions of eq.© searched. It must 
be noted that the presence of the Gamma function in the inverse Mellin transform turns 
into a very amenable numerical treatment, as pointed out in the following. 

The quite simple way we derived the latter relations entails many aspects about 
analytic functions, due to the complex nature of the inverse Mellin transform involved. 
A few remarks are therefore in order: i) Passing from eq. (112at to eq. (112bb it does not 
suffice a change of sign because of the presence of the integrand residue in y — 0. ii) 
It is not difficult to show that the integrand of eq.( |12at might have, at most, poles for 
<p{-k + iO), k — 0, 1,2... and, conversely, the integrand in eq.( |12bt might have, at most, 
poles for (p(k + iO), k = 1,2.... Hi) The integral in eq.(fT2l coincides with the sum of 
all the residues, and such sum can be proved to correspond to the r.h.s. of eq.©. Dis- 
tributions with divergent moments can therefore be represented as sum of the residues 
of the integrand eq.([T2l. Such an interesting theoretical aspect cannot be derived from 
classical applications of the inverse Mellin transform in probability, because in such 
works the starting point is the density and not the characteristic function. 

Of course, the fundamental strip depends on the integrability of the CF,following 
eqsdTOt. It has been proved in ll35l that the strictest fundamental strip associated 
with every absolute convergent CF, thus including a-stable random variables and other 
power-law distributions, is the interval < p < 1 if one uses eq.( 112a| i. Then, in the 
following, we always assume such a restriction on p. 

The density function can be restored from eq. dl2at by inverse Fourier transform, 
that explicitly reads 

p(x) = T- l {ip (m) ; x) = — j- F ( r ) (ll<p) (0) U -?e- iux du+ 

(2k) i Jp-ioo ( v Jo 
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+ {l 7 _ip)(0) J (-u)- y e-' u,x duj dy (13) 

Further simplification leads to 

i pp+ioo 

pW = t^2 rwrd-^k-ixr^XixF- 1 +<(ixro(-i^ 1 }dr (14) 

(2ny Jp-icc K ' 

Eq. (fl~4b represents the density function by the integral form of the Taylor approxi- 
mation and can be interpreted as the generalization of eq.©. We note that the density 
is expressed in integral form without the need of partitioning, and no restriction on the 
non-negativity of the density is needed. This fact is the consequence of having used 
the Mellin transform of the characteristic function. A further simplification is possible 
in virtue of the Hermitian property of the characteristic function, since ip(u) = <p(—u) 
where the overbar means conjugation. Evaluation of the CF for u > suffices to restore 
the PDF in the form 

p(x) = — ; Re _^ <( _ix r y>(i x )r-i dr ( i5) 

(2tt) 2 [J p-ioo sin (Try) J 

where the property nj sin(7ry) = F(l - y)F(y) has been used. By looking at eq.(TT3b 
one may state that the PDF in terms of complex moments remains meaningful and com- 
putationally useful. Moreover, since the integration is performed only on the imaginary 
axis, i.e. the real part of y remains fixed in the integrand, ((— iX) _y ) exists also for distri- 
butions for which (X 7 ) does not exist, for j integer and greater than a certain value. By 
comparing eqs.(fT2t and dl3T > a perfect duality between the representation of the PDF 
and the CF is evidenced. The integrals in eqs. dl2ab and (15[ remain independent of the 
choice of p within the strip < p < 1 , because of the holomorphy of the integrand in 
the fundamental strip. 

The appeal of such a representation appears evident when the discretization of the 
integral is performed. Indeed, it will be shown in the numerical section that the integral 
performed on the imaginary axis may be truncated at a certain value 77, because the 
integrand vanishes rapidly, and therefore the discretization of such an integral always 
produces very accurate results, in particular for distributions with heavy tails. 



3 Extension to Multivariate Random Vectors 

In this section, the concepts briefly outlined in the previous section will be derived 
in case of multivariate random vectors. Multivariate random vectors are fully char- 
acterized in probabilistic setting by the joint probability density function, or by its 
spectral counterpart, namely the joint characteristic function; alternatively also by joint 
moments or by joint cumulants of every order. Let (S, ZR, V) be a probability space, 
X : (S,M) -> (R d ,S d ), that is X T = (X u ...,X d ) denotes a d-dimensional random 
vector. As usual, by denoting x T = (x\, the joint distribution function of X is 

defined as F(x) = Pr(X\ < X\, ....,Xd < xj) and its derivative is called joint probability 
density function, namely p(x) = d d F (x) ldx\...dxd. Let the vector u r = (u\, mj), 
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with Uj € R; the joint characteristic function <p(u) is the multidimensional Fourier 
transform of the joint density, that is 

Xco y-»CO 
... I exp(iu r x) p(x)dxi...dx d (16) 
co U -co 

Analogously to the one dimensional case, we first define joint complex moments 

... (+ix 1 r 7l -(+ix d r 7d p(x)dx 1 ...dx d (17) 

It will be shown that they are multidimensional fractional derivatives of the joint 
characteristic function evaluated in zero. To this purpose, let x T -(x\, ...,x d ), y T = 
(j\, yd) be elements of R d and / (x) a function of d variables. As reported in ( ||38l , 
p. 462), the Liouville-type fractional integral of order y with Rey^ > 0, k — 1, d, is 
defined as 

^ w - fwfej I - 1 ff " -ervc * « , * t 

(18) 

where the choice of signs must be coherent between the left and the right hand sides. 
Let the function / (x) be defined on the whole space R d such that J \ f (x)| dx\ ...dx d <oo 

and indicate the multidimensional Fourier transform of / (x) as usual in the form 



(Tf) (u) = T{f (x) ; u}"= f e i(a ' «•♦ " : '/(f) df, ...d£, ( 1 9.) 

and the multidimensional Mellin transform as 

(M/)(y) = A4{/(x); r } f ... f 1 /(£)d£...d£ j (20) 

Jo Jo 

The multidimensional Fourier transform of the Liouville-type fractional integral of 
/ (x) is related to the Fourier transform of / (x) by the following expressions 

r {(ti...J) (*) ; "} = (+i"ir ri - w^V" (Ff) (u) (21a) 
^' ((£..±/) w ; u} = (±i Ul r yi ... (±m d y y " {r- l f) oo (21b) 



as reported in ( [3811 , p. 474). Now, we show that joint complex moments arise in a 
natural way once the fractional integral of the joint characteristic function is calculated 
in zero. Indeed, eq. (121bb specified for the joint characteristic function gives the relation 

{(/L ± ^)(u) ; x} = (±i Xl )-y ... (±ix rf r r " T- 1 (<p(u) ;x) (22) 

The multidimensional inverse Fourier transform of the joint characteristic function 
is the joint PDF and therefore, a multidimensional Fourier transform yields 

... (±ix l r n ...(±ix d )^p(x)exp(m T x)dx l ...dx d (23) 

co U -co 
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that, calculated in u r = (0, 0) gives the fundamental relation 

ll.„ ± <p(0, ...,0) = <(±iX,r y ' ...(±iX d )-™) (24) 

With the same reasoning, it is possible to relate the fractional derivatives of the 
joint characteristic function calculated in zero in terms of complex moments as 

Dl,. ± <p (0, 0) = <(±iXiF' ... (±iX d y) (25) 

Once these results are achieved, we derive a generalized integral form of the Tay- 
lor series. Such an expansion may be obtained by applying the commutative prop- 
erty of integro-differential operators and by the inverse Mellin transform also in the 
multidimensional case. To this aim, we interpret the Liouville-type multidimensional 
fractional integral as the multidimensional Mellin transform of the function / (x) 



CO CO 

r(y l )...r(y, / )(/L,/)<v l v,)- J... J er l A 



o o 

x/(xi + £ 1 ,...,x d + &)d£i...d&, (26) 

evaluated at fixed x T = (x\ , x d ). Writing eq. (126*11 for the joint characteristic function, 
evaluated in the point u T = (0, 0) 



CO CO 

r(n)...r( rd )(/L ± ^)(o,...,o) = j ... j u*-\. M y f l 



X 





Xip(+u\, +u c i)dui...du c i, (27) 

it is easy to recognize that the r.h.s is the Mellin transform, eq.(f20t. of the l.h.s.. Intro- 
ducing eq.d24l). the former assumes the relevant form 





x <p (+u\, +u d ) dui...du d . (28) 

By performing the inverse Mellin transform, the joint characteristic function is 
rewritten as 

pi+ico p^+ico 

(p(ui,...,u d ) = — — - f ... f r(yi)...r(y d )x 
{2mf J J 

pi-ico p^-ico 

x <(-sign(m) \X X )^ ... (-sign(u d ) iX d y y ") x 

x\ Ul \^ ...\u d \^ d ri ...dy d (29) 

The integrals are performed by choosing pi, ...,p d inside the fundamental multidi- 
mensional strip. Following [39 1, assume that ip(u\,..., u d ) is such that 

<p(u 1 ,...,u d ) = OQu 1 \-' 11 ,...,\u d \-i*) for nun(|Mi|,...,|M ( ,D->oo (30) 
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with qj > 0, j = 1, ...,d. Then, the multidimensional Mellin transform in eq.(t20b is 
meaningful under the conditions < pj < qj. Due to the absolute convergence of the 
joint characteristic function of a true multivariate random variable (see [35]) one can 
select the strictest fundamental strip < pj < 1 in which eq.j29t holds true. 

As in the one dimensional case, from the knowledge of the joint characteristic func- 
tion one can restore the joint PDF by Fourier inversion. The multidimensional inverse 
Fourier transform of relation (f28j) can be performed by subdivision of the domain R d 
in m = 2 d subsets in which the variables u T = (u\,.,.,ua) assume all the possible 
combinations of positive and negative values. By means of some algebra, the multidi- 
mensional inverse Fourier transform of eq.(|29b leads to the expression of the joint PDF 
of the random vector X , that reads 

pi+ioo p^+ioo 



(2m) d (2n) d J J sin (7171) sin (Try,/) 

Pi-ioo P,/-i°° 



2 

X 

h = \ />.;=! 

Vi- 1 I, ,^,+1 • \r<r-l 



^...^<((-l)Mx 1 )- r '...((-l)^iX d p>x 



^{(-\) k ^\x x f ...((-D^ +1 i^f d ri ...dy,, (31) 

where the relation nj sin (nyjj = F^l - ynTyyj) has been used. Eq.(l3TT> can be sim- 
plified recalling that the joint characteristic function has symmetry properties. Indeed, 
<p (u) = <p (-u) and the joint PDF can be restored by 2 d ~ l subsets, in the following form 

pi+ioo p ( /+ioo 



p(x u ...,x d ) =2Re 



> d (2n) d I "' / 



(2m) d (2n) d J J sin(7ryi)"""sin(7ry d ) 

pi-io 

2 

X ^ 

fe=l k d =\ 

x (i Xl r- 1 ((-lt +l ix 2 f^ ... ((-l)^ +1 ixaf d 7l ...dy d ) (32) 



: £ ... £ ((-iZO- {(-It iX 2 y 72 ... ((-If" ixf) x 



We note that in eq.([32b there occur d — 1 summations and the first variable is not 
indexed. Concluding, from eq.d29l> and d32b we can affirm that the function 

<(~iXi)-* ((-it- ix 2 y y2 ... ((-1/" ix d y y ") 

is able to represent both the joint PDF and the joint characteristic function. In the 
following section, it will be shown that discretization produces the desired results. 
That is, with a finite number of fractional moments, both the joint PDF and the joint 
characteristic function can be restored. 



4 Applications 

In this section it will be shown how previous results can be used to characterize proba- 
bility distributions. The integral forms in the monovariate case, eqs.(fT2l and (1141 1. and 
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in the multivariate case, eqs.d29i> and OTb . will be discretized in a suitable form. A 
remarkable effect of having applied the Mellin transform to the characteristic function 
rather than to the density is the presence of the gamma function in the integrand. In- 
deed, by its fast decay with increasing imaginary part of the argument, it makes the 
integral amenable to numerical calculations. 

4.1 Numerical treatment of the integral representation by complex 
moments for random variables 

The integrals in eqs.([T2b will be evaluated considering a partition of the truncated do- 
main [— 77, 77] in 2m (m e N) intervals of amplitude A77, such that t\ = mArj. In this way 
the integrals are approximated as a sum evaluated at yj = p + ikArj: 




(33a) 



2m+l 



(33b) 



fe=0 



and are valid for u > 0. For u < the characteristic function is restored by the prop- 
erty ip (u) = tf(-u). As mentioned, the truncation is feasible because of the presence of 
the gamma function that vanishes as the imaginary variable 77 increases. Consequently, 
according to eq.dTSl. the density can be evaluated by means of the same complex mo- 
ments in the form 



At? 
2n 



2m +1 



(ix) 



A=0 



sin {njk-m) 



(34) 



and, then, the calculation of 2m + 1 complex moments suffices to have a complete 
probabilistic description of the random variable. 



4.1.1 Example: Power-law tail distributions 

Firstly, we consider an c-stable random variable, with a = 1/2, whose characteristic 
function is known to be of the form [40 1 

(p(ll) = e _ °"" 2 l"l" 2 ( I_i sign(u))+i/j« 

where cr and fi are the scale and the location parameters, respectively. This distribution, 
also known as Levy-Smirnov distribution, has moments (X p ) < 00 only if p < 1 /2. 
By selecting A77 = .9, p = .85, and m = 60 the results plotted in the bi-logarithmic 
diagram in FigQ]show that in a very wide range the sum well approximates the integral. 
Complex moments have been calculated from eqs.dTOl 

((-iX)- y ) = (/V)(0) = -^—MW{+u)-y) 

y ' r (r) 
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by means of Mathematica. It has to be remarked that the results have been computed by 
applying eq.( ll2ab . and that p - .85 corresponds to complex moments /(— iX) -085 " 1 "^ 
which exist for this distribution. As further example, let us consider the 3/2-stable 
distribution, known as Holtsmark distribution, whose characteristic function is 

In this case, moments exist if p < 3/2 and, in order to apply eq.( 112a| >, the value of 
p = 0.35, i.e. moments of the type (i-iX)~° 35+in ) < oo are used. Other parameters in 
the discretization are A77 = .2 and m = 50 and the results have been reported in Figf2] 
We stress that the complex moments for both the examples have been calculated from 
the knowledge of the characteristic function without explicit use of the density function. 
This is quite interesting with respect to extensions of the method to random variables 
whose density is not known in explicit analytical form, such as a-stable densities. 



4.1.2 Statistics from raw data 

In this example it is shown that the proposed representation is also useful when dealing 
with raw data from experiments. Let us consider N s realizations of the random variable 
X and indicate by Xj, j = 1, ...,N S the j - th outcome of X. In order to evaluate the 
statistics, complex moments of order jk - p + i^t can be calculated from the raw data 
via 

for every k — 0, ...,2m + 1 following the partition of the previous example. Then, 
applying eq.d34ll the density can be computed. In Fig|3] a comparison between the 
histogram and the density calculated in such a way is shown, demonstrating excellent 
convergence. The parameters used in this example are N s - 10000, p = 0.5, A77 = 0.05, 
m = 50. 



4.2 Applications to multivariate random variables 

Joint densities and characteristic functions can be represented exactly by the integral 
forms in eqs.d29|i and d32| l in which joint complex moments appear. In this section, the 
integral is computed by a simple rectangle scheme, in order to present some applica- 
tions, without having to find the best numerical treatment of the integral involved. Yet, 
it has to be noted that more sophisticated numerical schemes (involving trapezoidal, 
Gaussian or other local/adaptive rules) of course may give better results. Hereinafter, 
the simplest integration scheme is adopted in order to show that with a finite number 
of quantities (the complex moments) both the joint PDF and the joint characteristic 
function are very well approximated. Then, indicating the integrand in eqs.d29l) and 
(f32l as 

G 9 (u,n,...,y d ) = — ^r(yi) ..Xte) ((-sign( Ml ) iX^ 1 " ... (-sign(w d ) \X d y 7i ) | Ml p 
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and 



G p (x,yi,...,y d ) - (2n[ f (2x) d S i„(^ yi) -.sin(^ yi ,) x 



X I - 2 (((-lf l iX i r\..((-l)^ixX 7i )x 

y- 1 t, ^jb+i • xt*- 1 



x((-D*' +1 ijt 1 f- 1 ...((-i)*' +1 i Jta ) 

we can more compactly express eqs.(|29i> and (l32t in the form 

Pi+ioo p rf +ioo 

¥>(u) = J" - J G ip (u,yi,...,y d )dyi...dy, 

pi-ioo Prf-i°° 

and 

pi+ioo Prf+io 



p(x)= J" ... J" G p (x,yi,...,y d )dyi...dy d 



pi-loo Pd-ico 



The former can now be approximated by the positions: yj = pj + b]j, rjj = rjAj, 
with — m < r j < m, m e N and Aj e R. Extremes of integration have been truncated up 
to [-mAy, mAyJ. The approximated dual series representing the joint statistics reads 

oo oo 

<p(u)=i d J] ... Yj G tp (vL,pi+inA u ...,Pd + ir d A d )Ai...A d 

n— -oo rd=-oo 



oo oo 

p(x)=i d ^ ... ^ G p (x,pi +inA u ...,Pd+ir d A d )Ai...A d 

ri=-oo ^=-00 

Moreover, the presence of the gamma function in the integrand of eqs.(|29il and (f32l > 
ensures that the integrand vanishes and, consequently, the approximated sums can be 
truncated in the form 

mi m c \ 

ip(u) = \ d ^ ... ^ G< p (u,pi+inA u ...,p d + ir d A d )A 1 ...A d 

n=-Mi r d =-m d 



p(x)=i d ^ ... ^ G p (x,pi +inA u ...,Pd+ir d A d ) Ai...A d 

n=-mi r d =-m d 

with m y e N. In this way, the accuracy of the results is of course affected by dis- 
cretization and a truncation error. By means of some applications we show that the 
truncated sums are nevertheless very good approximations, computing a finite number 
of complex moments. 
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4.2.1 Bivariate Cauchy random vector 



In this example we deal with the bivariate Cauchy distribution, that is an a-stable dis- 
tribution with a - I. Then, it has power-law tails and, as pointed out already in the 
introduction it has not integer moments of any order. The expression of its joint PDF 
is in general 

, , cdetff)-" 2 

P (X) = Jjr 

(l+ix-fifJl-Hx-^) 1 
while its joint characteristic function has the form 

(j- \l/2 -j- 
-u T,u ) + iu // 

and c = 7r~ 3 / 2 r(3/2). The matrix £ and the fi have the meaning of scale matrix and 
location vector, in some sense analogous to the covariance matrix and the mean vector. 
For the application of the fractional moments approach, we consider in this example 
the parameters 




then, as the location vector is not zero, we will expect a density with the maximum 
in the point x\ = 1 and X2 — 1 . Consequently, due to the condition fi + 0, the joint 
characteristic function of such random vector is a complex function as plainly descends 
from the definition. Then, eq.d29t and eq.d32li are discretized by choosing the steps 
A] = A2 = 0.3, and the integration is performed in the fundamental strip with p\ — 
P2 = 0.5. By using m = 20 fractional moments we show that very good numerical 
convergence is obtained. Indeed, in Fig|4]we first give the contour plots of the real and 
the imaginary part of the complex characteristic function where the dashed line is the 
approximation proposed. Without performing any further Fourier transform but just 
applying the discrete form of eq.(l32ll. with the same 20 fractional moments the joint 
PDF is calculated and plotted in Fig|5] showing a good numerical convergence in the 
whole domain. 



5 Conclusions 

In this paper, we provide new insights into the problem of the statistical representa- 
tion of random variables by means of the Mellin transform and fractional calculus. We 
define a class of complex moments (comprehending fractional moments) that can rep- 
resent every kind of joint PDF and joint characteristic function, including joint charac- 
teristic functions that are non-differentiable in zero, and densities with heavy-tails, like 
the bivariate Cauchy distribution or other power-law distributions. Comparison with 
classical applications of inverse Mellin transform shows that starting from the Mellin 
transform of the characteristic function and exploiting its Hermitian nature an easy and 
efficient representation of the statistics of random variables and vectors is achieved. 
Further, the complex moments introduced represent both the joint PDF and the joint 
CF. 
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The perspective presented in this paper might open, to the authors' opinion, a new 
touch in the moment problem, that is, under which conditions the knowledge of mo- 
ments suffices to reconstruct the underlying density, to hierarchy of moments and re- 
lated issues. 
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Figure 1 : Comparison between the exact (continuous) and the approximated (dotted) 
density of a 1/2-stable distribution with cr = 1 .0, fi = 2; power-law trend is plotted in 
dashed line. 
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Figure 2: (Left panel) comparison between the exact (continuous) and the approx- 
imated (dotted) characteristic function of a symmetric 3/2-stable distribution; (right 
panel) corresponding density in log-log diagram, contrasted with power-law trend 
(dashed line). 
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Figure 4: Comparison between the exact (continuous) and the approximated (dashed) 
real part of the joint characteristic function (left panel) and the imaginary part of the 
joint characteristic function(right panel) of bivariate Cauchy distribution 
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Figure 5: Comparison between the exact (continuous) and the approximated (dashed) 
joint PDF of bivariate Cauchy distribution with no-zero location 
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